Index
Murray Logan
31 August, 2019
Preparations
library(caret) # load this first as it loads plyr
library(gbm)
library(car)
library(grid)
library(gridExtra)
library(randomForest)
library(xtable)
library(kableExtra)
library(abt)
library(tidyverse)
source('FISH_functions.R')Raw fish data
fish = read.csv('data/Fish benthic physical sitelevel 2019.csv', strip.white=TRUE)
fish %>% glimpse## Observations: 619
## Variables: 300
## $ YEAR <int> 2007, 2007, 2007, 2007, 2007, 2007, …
## $ REGION <fct> Keppel, Keppel, Keppel, Keppel, Kepp…
## $ EXPOSURE <fct> Semi-Exposed, Exposed, Semi-Exposed,…
## $ NTR <fct> Fished, Fished, Fished, Fished, Fish…
## $ NTR.Pooled <fct> Fished, Fished, Fished, Fished, Fish…
## $ SITE <fct> GK1, GK2, GK3, GK4, GK7, GK8, H3, MI…
## $ lut.arge <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ lut.carp <dbl> 13.3320, 11.3322, 4.6662, 8.6658, 16…
## $ lut.flma <dbl> 0.0000, 3.3330, 0.0000, 0.0000, 7.33…
## $ lut.fulv <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ lut.kas <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ lut.lemn <dbl> 0.0000, 0.6666, 0.0000, 0.0000, 0.00…
## $ lut.lutj <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ lut.mono <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ lut.quin <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ lut.rivu <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ lut.russ <dbl> 67.3266, 4.6662, 0.0000, 5.3328, 116…
## $ lut.seba <dbl> 1.9998, 0.6666, 9.9990, 0.0000, 0.00…
## $ lut.vitt <dbl> 377.9622, 10.6656, 16.6650, 9.3324, …
## $ mac.macu <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ pagrus <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ sym.nema <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ par.barb <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ par.boides <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ par.bifa <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ par.indi <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ par.multi <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ par.cili <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ upe.trag <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ platax <dbl> 0.0000, 0.6666, 0.0000, 0.6666, 0.66…
## $ Platycephalus.spp. <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ let.atki <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ let.hara <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ let.lati <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ let.lent <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ let.mini <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ let.nebu <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ let.oliv <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ let.obso <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ let.orna <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ gym.spp <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ mon.gran <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ sco.bili <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ sco.marg <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ sco.mono <dbl> 4.6662, 1.9998, 3.9996, 5.3328, 3.99…
## $ sco.line <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ any.leuc <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ ath.roga <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cep.argu <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cep.boen <dbl> 21.9978, 5.9994, 7.3326, 8.6658, 3.3…
## $ cep.cyan <dbl> 1.3332, 0.0000, 1.9998, 1.3332, 0.00…
## $ cep.micr <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cro.alti <dbl> 0.6666, 0.0000, 0.6666, 0.6666, 0.00…
## $ epi.caer <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ epi.coio <dbl> 0.0000, 0.0000, 0.0000, 1.3332, 0.00…
## $ epi.cora <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ epi.fasc <dbl> 0.0000, 0.0000, 0.6666, 0.0000, 0.00…
## $ epi.fusc <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ epi.howl <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ epi.lanc <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ epi.hexa <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ epi.merr <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ epi.ongu <dbl> 0.6666, 0.0000, 0.6666, 0.0000, 0.00…
## $ epi.quoy <dbl> 6.6660, 2.6664, 12.6654, 1.9998, 1.9…
## $ epi.sexf <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ pms.laev <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ pms.leop <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.66…
## $ pms.macu <dbl> 5.9994, 3.9996, 17.3316, 7.3326, 9.9…
## $ pte.ante <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ pte.voli <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ dip.bifa <dbl> 1.9998, 1.9998, 7.9992, 1.9998, 1.99…
## $ dia.pict <dbl> 23.3310, 3.3330, 4.6662, 1.9998, 9.9…
## $ ple.chae <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ ple.flav <dbl> 1.9998, 1.3332, 2.6664, 0.6666, 0.00…
## $ ple.gibb <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ ple.less <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ ple.line <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ ple.picu <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ ple.unic <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ bodianus <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ oxy.diag <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ che.chlo <dbl> 1.3332, 0.6666, 5.9994, 0.6666, 1.33…
## $ che.fasc <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ che.tril <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ che.undu <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cho.anch <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cho.ceph <dbl> 0.6666, 0.0000, 0.6666, 1.3332, 0.66…
## $ cho.cya <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ cho.fasc <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ cho.grap <dbl> 0.0000, 0.0000, 0.6666, 0.0000, 1.33…
## $ cho.scho <dbl> 1.3332, 0.0000, 0.6666, 0.6666, 0.00…
## $ cho.vitt <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ epb.insi <dbl> 0.6666, 0.0000, 0.0000, 0.0000, 1.33…
## $ gom.vari <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ hem.melt <dbl> 17.3316, 9.3324, 25.9974, 13.9986, 2…
## $ hem.fasc <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ psa.wai <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ aca.bloc <dbl> 0.0000, 0.0000, 0.6666, 0.0000, 4.66…
## $ aca.duss <dbl> 2.6664, 0.0000, 0.0000, 0.0000, 9.99…
## $ aca.gram <dbl> 0.0000, 0.0000, 0.6666, 1.3332, 4.66…
## $ aca.line <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ aca.nans <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ aca.ncus <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ aca.nuda <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ aca.xant <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ aca.thom <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ aca.spp <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cte.bino <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cte.stri <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ nas.annu <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ nas.brach <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ nas.brevi <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ nas.litu <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ nas.tong <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ nas.unic <dbl> 0.0000, 1.9998, 0.0000, 0.0000, 0.00…
## $ pri.spp <dbl> 0.0000, 0.0000, 0.0000, 3.9996, 18.6…
## $ kyp.spp <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ zeb.scop <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ zeb.veli <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ zan.corn <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ bol.muri <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cal.caro <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cet.bico <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ chl.blee <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ chl.micr <dbl> 0.6666, 0.0000, 0.0000, 0.0000, 0.00…
## $ chs.sord <dbl> 10.6656, 2.6664, 0.0000, 0.0000, 0.0…
## $ hip.long <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ sca.alti <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ sca.cham <dbl> 6.6660, 0.6666, 0.0000, 0.0000, 0.00…
## $ sca.dimi <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ sca.flav <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ sca.fors <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ sca.fren <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ sca.ghob <dbl> 7.9992, 3.9996, 1.9998, 0.6666, 13.3…
## $ sca.glob <dbl> 0.6666, 7.9992, 0.0000, 0.0000, 0.66…
## $ sca.nigr <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ sca.ovic <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ sca.psit <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ sca.rivu <dbl> 37.3296, 25.3308, 19.3314, 7.9992, 4…
## $ sca.rubr <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ sca.schl <dbl> 0.6666, 0.0000, 0.0000, 0.0000, 0.00…
## $ sca.spin <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ sca.tric <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ sca.spp <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ sig.arge <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ sig.cana <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ sig.cora <dbl> 0.0000, 1.3332, 0.0000, 0.0000, 0.00…
## $ sig.doli <dbl> 5.3328, 3.9996, 1.3332, 2.6664, 6.66…
## $ sig.fusc <dbl> 0.0000, 0.0000, 11.9988, 6.6660, 0.0…
## $ sig.javu <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ sig.line <dbl> 0.0000, 0.0000, 0.0000, 1.3332, 0.00…
## $ sig.puel <dbl> 0.0000, 0.0000, 0.6666, 0.0000, 0.00…
## $ sig.ptus <dbl> 0.0000, 0.0000, 0.0000, 1.9998, 0.00…
## $ sig..ptiss <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ sig.stell <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ sig.vulp <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ mon.arge <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cha.auri <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ cha.afas <dbl> 63.3270, 89.9910, 39.3294, 57.3276, …
## $ cha.baro <dbl> 0.0000, 1.3332, 0.0000, 0.0000, 0.00…
## $ cha.benn <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cha.ephi <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ cha.citr <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cha.flav <dbl> 0.0000, 0.0000, 0.0000, 0.6666, 0.00…
## $ cha.line <dbl> 5.9994, 5.3328, 0.0000, 1.3332, 5.33…
## $ cha.lunu <dbl> 1.3332, 0.0000, 0.6666, 1.3332, 0.00…
## $ cha.ttus <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cha.mela <dbl> 4.6662, 1.3332, 1.9998, 0.0000, 0.66…
## $ cha.ocel <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cha.orna <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cha.pleb <dbl> 0.0000, 0.0000, 0.6666, 0.6666, 0.00…
## $ cha.rain <dbl> 3.3330, 13.9986, 1.3332, 3.3330, 3.9…
## $ cha.raff <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cha.seme <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cha.spec <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ cha.tlis <dbl> 1.3332, 2.6664, 0.0000, 0.0000, 0.00…
## $ cha.ulie <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ cha.vaga <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ chm.rost <dbl> 7.3326, 3.3330, 1.9998, 1.3332, 6.66…
## $ cor.alti <dbl> 3.9996, 1.9998, 0.0000, 0.6666, 1.99…
## $ cor.chry <dbl> 0.0000, 1.3332, 0.0000, 0.6666, 1.99…
## $ hen.acum <dbl> 0.0000, 1.3332, 0.0000, 0.0000, 0.00…
## $ hen.mono <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ hen.vari <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ mct.stri <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.00…
## $ prc.ocel <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cen.bico <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cen.bisp <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cen.nox <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cen.tibi <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ cen.vrol <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ chd.doub <dbl> 5.9994, 3.9996, 0.6666, 1.3332, 7.99…
## $ chd.mere <dbl> 5.3328, 5.3328, 1.9998, 2.6664, 4.66…
## $ pmc.impe <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ pmc.semi <dbl> 0.0000, 3.3330, 0.0000, 0.0000, 0.66…
## $ pmc.sext <dbl> 1.3332, 1.3332, 0.0000, 0.6666, 0.66…
## $ pmc.xant <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ pyg.diac <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ abu.spp <int> 406, 48, 0, 0, 238, 0, 30, 112, 38, …
## $ acn.poly <int> 90, 48, 74, 58, 124, 22, 152, 78, 62…
## $ amb.aure <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ amb.cura <int> 0, 0, 4, 2, 0, 0, 0, 0, 0, 2, 0, 0, …
## $ amb.leuc <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ amp.spp <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ amp.peri <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ pre.spp <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ chl.labi <int> 8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ chr.ambo <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ chr.alis <int> 60, 0, 12, 0, 0, 0, 138, 34, 16, 64,…
## $ chr.apes <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ chr.marg <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ chr.niti <int> 8920, 9680, 7340, 15840, 12210, 240,…
## $ chr.retr <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ chr.tern <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ chr.webe <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ chr.xant <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ chy.roll <int> 0, 0, 2, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ chy.talb <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ chy.flav <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ das.spp <int> 12, 0, 0, 0, 0, 0, 0, 0, 0, 24, 4, 0…
## $ das.reti <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ das.trim <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ dis.spp <int> 0, 0, 8, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ hgy.plag <int> 12, 0, 0, 0, 0, 0, 0, 0, 0, 0, 6, 0,…
## $ neg.mela <int> 14, 6, 2, 4, 12, 0, 4, 0, 8, 0, 0, 2…
## $ parma <int> 10, 14, 0, 0, 4, 2, 2, 2, 4, 0, 4, 0…
## $ neg.nigr <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ pgy.dick <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ pgy.lacr <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ pom.adel <int> 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, …
## $ pom.ambo <int> 0, 0, 8, 0, 0, 0, 18, 6, 4, 32, 0, 0…
## $ pom.aust <int> 58, 140, 14, 54, 28, 0, 216, 156, 21…
## $ pom.bank <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ pom.brac <int> 84, 0, 4, 6, 52, 0, 90, 14, 10, 0, 6…
## $ pom.chry <int> 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, …
## $ pom.coel <int> 0, 0, 2, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ pom.gram <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ pom.lepi <int> 8, 0, 6, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ pom.molu <int> 602, 294, 426, 634, 542, 190, 564, 1…
## $ pom.naga <int> 6, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ pom.phil <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ pom.reid <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ pom.trip <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ pom.vaiu <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ pom.ward <int> 258, 224, 206, 338, 200, 456, 136, 1…
## $ ste.apic <int> 2, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, …
## $ ste.fasc <int> 8, 8, 4, 10, 8, 12, 0, 4, 8, 0, 4, 6…
## $ ste.livi <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ ste.nigr <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ oxy.long <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ Anampses <int> 0, 6, 0, 2, 8, 0, 0, 0, 0, 0, 6, 2, …
## $ Coris <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ Halichoeres <int> 14, 10, 20, 30, 8, 8, 26, 6, 10, 48,…
## $ Labroides <int> 18, 10, 10, 12, 16, 4, 10, 10, 4, 2,…
## $ Labropsis <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ Psuedolabrus.guentheri <int> 20, 58, 20, 24, 26, 104, 16, 22, 20,…
## $ Stethojoulis <int> 12, 4, 16, 28, 2, 30, 16, 2, 0, 84, …
## $ Thalasoma <int> 24, 8, 12, 12, 24, 2, 30, 10, 16, 18…
## $ Macropharyngodon <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ Labrichthys <int> 8, 0, 0, 6, 0, 2, 0, 2, 4, 4, 4, 0, …
## $ Total.fish.density <dbl> 11381.9272, 10805.9752, 8401.3122, 1…
## $ Total.fish.species.richness <dbl> 35.2, 29.4, 27.6, 26.4, 29.2, 17.8, …
## $ Prey.density <int> 10654, 10558, 8190, 17060, 13502, 10…
## $ Prey.biomass <dbl> 44.760041, 30.942239, 23.789262, 46.…
## $ Plectropomus.total.density <dbl> 5.9994, 3.9996, 17.3316, 7.3326, 10.…
## $ Plectropomus.total.biomass <dbl> 3.980542, 5.371889, 9.546584, 8.9263…
## $ Plectropomus.legal.density <dbl> 1.9998, 2.6664, 4.6662, 3.9996, 1.99…
## $ Plectropomus.legal.biomass <dbl> 2.4194697, 4.5827655, 5.0227904, 7.1…
## $ BE <dbl> 121.9966, 119.9976, 114.6630, 129.33…
## $ GRAZ <dbl> 72.6594, 47.9952, 36.6630, 26.6640, …
## $ COR <dbl> 103.3246, 119.3214, 45.9954, 71.9934…
## $ OM <dbl> 499.9996, 99.9996, 74.0000, 59.9998,…
## $ PL <dbl> 9000, 9680, 7366, 15842, 12210, 240,…
## $ CA <dbl> 521.9478, 45.9954, 65.9934, 43.9956,…
## $ PI <dbl> 7.9992, 6.6660, 25.9974, 10.6656, 12…
## $ FA <int> 290, 246, 218, 348, 212, 474, 138, 1…
## $ slope <dbl> 2.60, 2.64, 2.00, 2.40, 1.96, 2.20, …
## $ rugosity <dbl> 2.68, 3.08, 3.00, 3.00, 3.00, 2.68, …
## $ SCI <dbl> 7.040, 8.088, 6.000, 7.200, 5.880, 5…
## $ LCC_. <dbl> 60.0, 78.8, 77.6, 77.6, 57.2, 13.6, …
## $ LHC_. <dbl> 56.4, 68.8, 77.6, 75.6, 55.6, 13.6, …
## $ SC_. <dbl> 3.6, 10.0, 0.0, 2.0, 1.6, 0.0, 0.0, …
## $ MAC_. <dbl> 2.4, 0.4, 9.6, 4.4, 3.6, 4.0, 29.2, …
## $ SPO_. <dbl> 1.6, 0.0, 0.0, 0.4, 0.0, 0.0, 0.0, 0…
## $ Turf_. <dbl> 26.0, 17.6, 4.4, 11.6, 28.0, 71.2, 6…
## $ Unconsolidated_. <dbl> 9.6, 3.2, 8.4, 6.0, 11.2, 11.2, 7.2,…
## $ Benthic.richness <dbl> 6.6, 5.8, 2.2, 6.0, 3.8, 2.0, 3.4, 5…
## $ Coral_Morph.Diversity <dbl> 5.4, 5.6, 1.2, 5.0, 3.4, 1.4, 2.4, 4…
## $ ChlA <dbl> 0.45221510, 0.45221510, 0.41940367, …
## $ kd490 <dbl> 0.06519040, 0.06519040, 0.06261333, …
## $ SSTmean <fct> 23.80386083, 23.80357108, 23.7912496…
## $ SSTanom <dbl> -0.3130472, -0.3140394, -0.3159858, …
## $ wave.exposure.index <dbl> 0.007256172, 0.160425811, 0.00099007…
## $ Corrected.depth <dbl> 4.03, 6.39, 5.12, 5.47, 3.73, 3.02, …
## $ maxDHW <dbl> 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, …
## $ deltaLHC <dbl> 14.0, 7.6, 6.4, 16.8, 8.0, 2.8, 4.8,…
## $ deltaMA <dbl> -6.8, -1.2, -3.6, 0.4, -6.4, -41.2, …
## $ Cyclone <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ LT.Fprimary <dbl> 0.087, 0.087, 0.087, 0.031, 0.025, 0…
## $ Exposure.to.primary.weeks <int> 2, 2, 2, 2, 3, 2, 3, 1, 1, 0, 3, 2, …var.lookup = rbind(
data.frame(pretty.name='Total density', Field.name='Total.fish.density', Abbreviation='TFD', Family='gaussian', Type='Response', Transform='log', Groupby=''),
data.frame(pretty.name='Total species richness', Field.name='Total.fish.species.richness', Abbreviation='TFSR', Family='gaussian', Type='Response', Transform='I', Groupby=''),
data.frame(pretty.name='Benthic invertivores', Field.name='BE', Abbreviation='BE', Family='gaussian', Type='Response', Transform='log', Groupby=''),
data.frame(pretty.name='Grazers', Field.name='GRAZ', Abbreviation='GRAZ', Family='gaussian', Type='Response', Transform='log', Groupby=''),
data.frame(pretty.name='Corallivores', Field.name='COR', Abbreviation='COR', Family='gaussian', Type='Response', Transform='log', Groupby=''),
data.frame(pretty.name='Omnivores', Field.name='OM', Abbreviation='OM', Family='gaussian', Type='Response', Transform='log', Groupby=''),
data.frame(pretty.name='Planktivores', Field.name='PL', Abbreviation='PL', Family='gaussian', Type='Response', Transform='log', Groupby=''),
data.frame(pretty.name='Carnivores', Field.name='CA', Abbreviation='CA', Family='gaussian', Type='Response', Transform='log', Groupby=''),
data.frame(pretty.name='Piscivores', Field.name='PI', Abbreviation='PI', Family='gaussian', Type='Response', Transform='log', Groupby=''),
data.frame(pretty.name='Farmers', Field.name='FA', Abbreviation='FA', Family='gaussian', Type='Response', Transform='log', Groupby=''),
data.frame(pretty.name='Plectropomus total density', Field.name='Plectropomus.total.density', Abbreviation='PTD', Family='gaussian', Type='Response', Transform='log', Groupby=''),
data.frame(pretty.name='Plectropomus total biomass', Field.name='Plectropomus.total.biomass', Abbreviation='PTB', Family='gaussian', Type='Response', Transform='log', Groupby=''),
data.frame(pretty.name='Plectropomus legal density', Field.name='Plectropomus.legal.density', Abbreviation='PLD', Family='gaussian', Type='Response', Transform='log', Groupby=''),
data.frame(pretty.name='Plectropomus legal biomass', Field.name='Plectropomus.legal.biomass', Abbreviation='PLB', Family='gaussian', Type='Response', Transform='log', Groupby=''),
data.frame(pretty.name='Region', Field.name='REGION', Abbreviation='REGION', Family=NA, Type='Predictor', Transform='I', Groupby=''),
data.frame(pretty.name='NTR Pooled', Field.name='NTR.Pooled', Abbreviation='NTR.Pooled', Family=NA, Type='Predictor', Transform='I', Groupby='Region'),
data.frame(pretty.name='LHC % (live hard coral cover)', Field.name='LHC_.', Abbreviation='LHC', Family=NA, Type='Predictor', Transform='I', Groupby='Region'),
data.frame(pretty.name='SC % (soft coral)', Field.name='SC_.', Abbreviation='SC', Family=NA, Type='Predictor', Transform='I', Groupby='Region'),
data.frame(pretty.name='MA % (macroalgal cover)', Field.name='MAC_.', Abbreviation='MA', Family=NA, Type='Predictor', Transform='I', Groupby='Region'),
data.frame(pretty.name='Turf %', Field.name='Turf_.', Abbreviation='TURF', Family=NA, Type='Predictor', Transform='I', Groupby='Region'),
data.frame(pretty.name='Unconsolidated %', Field.name='Unconsolidated_.', Abbreviation='UC', Family=NA, Type='Predictor', Transform='I', Groupby='Region'),
data.frame(pretty.name='Benthic richness', Field.name='Benthic.richness', Abbreviation='BR', Family=NA, Type='Predictor', Transform='I', Groupby='Region'),
data.frame(pretty.name='Coral morphological diversity', Field.name='Coral_Morph.Diversity', Abbreviation='CMD', Family=NA, Type='Predictor', Transform='I', Groupby='Region'),
data.frame(pretty.name='Slope', Field.name='slope', Abbreviation='SLOPE', Family=NA, Type='Predictor', Transform='I', Groupby='Region'),
data.frame(pretty.name='Rugosity', Field.name='rugosity', Abbreviation='RUG', Family=NA, Type='Predictor', Transform='I', Groupby='Region'),
data.frame(pretty.name='SCI (structural complexity)', Field.name='SCI', Abbreviation='SCI', Family=NA, Type='Predictor', Transform='I', Groupby='Region'),
data.frame(pretty.name='Chla', Field.name='ChlA', Abbreviation='CHL', Family=NA, Type='Predictor', Transform='I', Groupby='Region'),
data.frame(pretty.name='Kd490', Field.name='kd490', Abbreviation='KD490', Family=NA, Type='Predictor', Transform='I', Groupby='Region'),
data.frame(pretty.name='SST mean', Field.name='SSTmean', Abbreviation='SSTMEAN', Family=NA, Type='Predictor', Transform='I', Groupby='Region'),
data.frame(pretty.name='SST anom', Field.name='SSTanom', Abbreviation='SSTANOM', Family=NA, Type='Predictor', Transform='I', Groupby='Region'),
data.frame(pretty.name='Wave exposure', Field.name='wave.exposure.index', Abbreviation='WAVE', Family=NA, Type='Predictor', Transform='I', Groupby='Region'),
data.frame(pretty.name='Corrected depth', Field.name='Corrected.depth', Abbreviation='DEPTH', Family=NA, Type='Predictor', Transform='I', Groupby='Region'),
data.frame(pretty.name='Max DHW', Field.name='maxDHW', Abbreviation='DHW', Family=NA, Type='Predictor', Transform='I', Groupby='Region'),
data.frame(pretty.name='Cyclone', Field.name='Cyclone', Abbreviation='CYCLONE', Family=NA, Type='Predictor', Transform='I', Groupby='Region'),
data.frame(pretty.name='Exposure to primary weeks', Field.name='Exposure.to.primary.weeks', Abbreviation='EXP', Family=NA, Type='Predictor', Transform='I', Groupby='Region'),
data.frame(pretty.name='Prey density', Field.name='Prey.density', Abbreviation='PREY.DENSITY', Family=NA, Type='Predictor', Transform='I', Groupby='Region'),
data.frame(pretty.name='Prey biomass', Field.name='Prey.biomass', Abbreviation='PREY.BIOMASS', Family=NA, Type='Predictor', Transform='I', Groupby='Region')
)
save(var.lookup, file='data/var.lookup.RData')
names = with(var.lookup, setNames(as.character(Field.name), Abbreviation))
## exclude those whose names equal their values otherwise there will be duplicate fields created in the fish data
names = names[names(names)!=names]## Create duplicates of the fields that are not already abbreviated
fish = fish %>%
mutate(SSTmean=ifelse(SSTmean=='#N/A',NA,SSTmean)) %>%
mutate_at(as.character(var.lookup$Field.name), list(A=~I)) %>%
rename(!!! gsub('(.*)','\\1_A',names)) %>%
mutate(REGION=factor(REGION, levels=c('Palm','Magnetic','Whitsunday','Keppel')))
save(fish, file='data/fish.RData')Response Variables:
| Variable | Data field | Abbreviation | Distribution |
|---|---|---|---|
| Total density | Total.fish.density |
TFD | Gaussian (log) |
| Total species richness | Total.fish.species.richness |
TFSR | Gaussian |
| Benthic invertivores | BE |
BE | Gaussian (log) |
| Grazers | GRAZ |
GRAZ | Gaussian (log) |
| Corallivores | COR |
COR | Gaussian (log) |
| Omnivores | OM |
OM | Gaussian (log) |
| Planktivores | PL |
PL | Gaussian (log) |
| Carnivores | CA |
CA | Gaussian (log) |
| Piscivores | PI |
PI | Gaussian (log) |
| Farmers | FA |
FA | Gaussian (log) |
| Plectropomus total density | Plectropomus.total.density |
PTD | Gaussian (log) |
| Plectropomus total biomass | Plectropomus.total.biomass |
PTB | Gaussian (log) |
| Plectropomus legal density | Plectropomus.legal.density |
PLD | Gaussian (log) |
| Plectropomus legal biomass | Plectropomus.legal.biomass |
PLB | Gaussian (log) |
Exploratory data analyses
EDA_histograms = function(var='', group='',dat=NULL, var.lookup) {
pretty.name = as.character(var.lookup$pretty.name[var.lookup$Abbreviation==var])
if (is.numeric(dat[,var])) {
grid.arrange(
ggplot(dat) + geom_histogram(aes(x=!!sym(var))) +
scale_x_continuous(pretty.name),
ggplot(dat) + geom_histogram(aes(x=!!sym(var))) +
scale_x_log10(pretty.name),
nrow=1
)
} else {
grid.arrange(
ggplot(dat) + geom_bar(aes(x=!!sym(var), fill=!!sym(group))) +
scale_x_discrete(pretty.name),
nrow=1
)
}
}EDA_density = function(var='', group='', dat=NULL, var.lookup) {
pretty.name = as.character(var.lookup$pretty.name[var.lookup$Abbreviation==var])
if (is.numeric(dat[,var])) {
grid.arrange(
ggplot(dat) + geom_density(aes(x=!!sym(var), fill=!!sym(group)), alpha=0.5) +
scale_x_continuous(pretty.name),
ggplot(dat) + geom_density(aes(x=!!sym(var), fill=!!sym(group)), alpha=0.5) +
scale_x_log10(pretty.name),
nrow=1
)
} else {
}
}Total density
EDA_histograms(var='TFD', dat=fish, var.lookup=var.lookup)
EDA_density(var='TFD', group='REGION', dat=fish, var.lookup=var.lookup)
Total species richness
EDA_histograms(var='TFSR', dat=fish, var.lookup=var.lookup)
EDA_density(var='TFSR', group='REGION', dat=fish, var.lookup=var.lookup)
Benthic invertivores
EDA_histograms(var='BE', dat=fish, var.lookup=var.lookup)
EDA_density(var='BE', group='REGION', dat=fish, var.lookup=var.lookup)
Grazers
EDA_histograms(var='GRAZ', dat=fish, var.lookup=var.lookup)
EDA_density(var='GRAZ', group='REGION', dat=fish, var.lookup=var.lookup)
Corallivores
EDA_histograms(var='COR', dat=fish, var.lookup=var.lookup)
EDA_density(var='COR', group='REGION', dat=fish, var.lookup=var.lookup)
Omnivores
EDA_histograms(var='OM', dat=fish, var.lookup=var.lookup)
EDA_density(var='OM', group='REGION', dat=fish, var.lookup=var.lookup)
Planktivores
EDA_histograms(var='PL', dat=fish, var.lookup=var.lookup)
EDA_density(var='PL', group='REGION', dat=fish, var.lookup=var.lookup)
Carnivores
EDA_histograms(var='CA', dat=fish, var.lookup=var.lookup)
EDA_density(var='CA', group='REGION', dat=fish, var.lookup=var.lookup)
Piscivores
EDA_histograms(var='PI', dat=fish, var.lookup=var.lookup)
EDA_density(var='PI', group='REGION', dat=fish, var.lookup=var.lookup)
Farmers
EDA_histograms(var='FA', dat=fish, var.lookup=var.lookup)
EDA_density(var='FA', group='REGION', dat=fish, var.lookup=var.lookup)
Plectropomus total density
EDA_histograms(var='PTD', dat=fish, var.lookup=var.lookup)
EDA_density(var='PTD', group='REGION', dat=fish, var.lookup=var.lookup)
Plectropomus total biomass
EDA_histograms(var='PTB', dat=fish, var.lookup=var.lookup)
EDA_density(var='PTB', group='REGION', dat=fish, var.lookup=var.lookup)
Plectropomus legal density
EDA_histograms(var='PLD', dat=fish, var.lookup=var.lookup)
EDA_density(var='PLD', group='REGION', dat=fish, var.lookup=var.lookup)
Plectropomus legal biomass
EDA_histograms(var='PLB', dat=fish, var.lookup=var.lookup)
EDA_density(var='PLB', group='REGION', dat=fish, var.lookup=var.lookup)
Predictor Variables:
| Variable | Data field | Abbreviation |
|---|---|---|
| Region | REGION |
REGION |
| NTR Pooled | NTR.Pooled |
NTR.Pooled |
| LHC % (live hard coral cover) | LHC_. |
LHC |
| SC % (soft coral) | SC_. |
SC |
| MA % (macroalgal cover) | MAC_. |
MA |
| Turf % | Turf_. |
TURF |
| Unconsolidated % | Unconsolidated_. |
UC |
| Benthic richness | Benthic.richness |
BR |
| Coral morphological diversity | Coral_Morph.Diversity |
CMD |
| Slope | slope |
SLOPE |
| Rugosity | rugosity |
RUG |
| SCI (structural complexity) | SCI |
SCI |
| Chla | ChlA |
CHL |
| Kd490 | kd490 |
KD490 |
| SST mean | SSTmean |
SSTMEAN |
| SST anom | SSTanom |
SSTANOM |
| Wave exposure | wave.exposure.index |
WAVE |
| Corrected depth | Corrected.depth |
DEPTH |
| Max DHW | maxDHW |
DHW |
| Cyclone | Cyclone |
CYCLONE |
| Exposure to primary weeks | Exposure.to.primary.weeks |
EXP |
FOR Plectropomus response variables ONLY, add:
| Variable | Data field | Abbreviation |
|---|---|---|
| Prey density (for Plectropomus density and legal density) | Prey.density |
PREY.DENSITY |
| Prey biomass (for Plectropomus biomass and legal biomass) | Prey.biomass |
PREY.BIOMASS |
Species-specific analyses: list of species in the attached file. There are two lists there: the single-column list is for the initial analysis between regions, and the table is for the species lists to use for each region individually.
Raw selections data
selections = read.csv('data/Species selection.csv', strip.white=TRUE)
selections %>% glimpse## Observations: 264
## Variables: 50
## $ Species <fct> aca.bloc, aca.duss, aca.gram, aca.line, …
## $ Palms <dbl> 2517.48, 3359.97, 1045.62, 0.00, 13.32, …
## $ Magnetic <dbl> 579.42, 679.32, 126.54, 0.00, 0.00, 0.00…
## $ Whitsundays <dbl> 1728.27, 965.70, 3716.28, 13.32, 76.59, …
## $ Keppels <dbl> 789.21, 346.32, 1312.02, 0.00, 23.31, 0.…
## $ Palms...total <dbl> 0.150705201, 0.201139613, 0.062594488, 0…
## $ Magnetic...total <dbl> 0.94095997, 1.10319444, 0.20549700, 0.00…
## $ Whitsundays...total <dbl> 0.053913389, 0.030125015, 0.115929368, 0…
## $ Keppels...total <dbl> 0.0145102620, 0.0063673720, 0.0241225460…
## $ Palms...FG <dbl> 2.816481633, 3.759034349, 1.169808509, 0…
## $ Magnetic...FG <dbl> 12.3404255, 14.4680851, 2.6950355, 0.000…
## $ Whitsundays...FG <dbl> 1.364317447, 0.762335375, 2.933676822, 0…
## $ Keppels...FG <dbl> 1.697099893, 0.744718940, 2.821339062, 0…
## $ X <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
## $ X.1 <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
## $ Species.1 <fct> pom.molu, pom.brac, pom.ward, acn.poly, …
## $ Palms...total.1 <dbl> 23.9783308, 16.8312258, 7.5841087, 5.090…
## $ Species.2 <fct> pom.ward, pom.molu, Halichoeres, lut.flm…
## $ Magnetic...total.1 <dbl> 38.3094226, 11.1891446, 7.4864959, 6.451…
## $ Species.3 <fct> pom.molu, pom.brac, chr.niti, acn.poly, …
## $ Whitsundays...total.1 <dbl> 21.6939414, 15.6405509, 15.1548448, 8.69…
## $ Species.4 <fct> chr.niti, pom.ward, pom.molu, pom.aust, …
## $ Keppels...total.1 <dbl> 85.50374438, 4.50911903, 2.86450864, 1.4…
## $ X.2 <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
## $ X.3 <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
## $ Species.5 <fct> sca.rivu, sig.doli, sca.flav, chs.sord, …
## $ Palms...FG.1 <dbl> 37.6387750, 19.0187020, 7.9837568, 4.157…
## $ Species.6 <fct> sca.rivu, sca.ghob, aca.duss, aca.bloc, …
## $ Magnetic...FG.1 <dbl> 39.92907801, 17.23404255, 14.46808511, 1…
## $ Species.7 <fct> sca.rivu, sig.doli, sca.flav, sca.nigr, …
## $ Whitsundays...FG.1 <dbl> 31.6789779, 13.0411924, 12.1106175, 10.3…
## $ Species.8 <fct> sca.rivu, sca.ghob, sig.arge, sig.fusc, …
## $ Keppels...FG.1 <dbl> 47.26817043, 15.33834586, 11.57178661, 6…
## $ X.4 <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
## $ X.5 <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
## $ X.6 <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
## $ X5.most.abundant <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
## $ FG <fct> Grazers, , , , , Planktivores, , , , , C…
## $ Palms.1 <fct> sca.rivu, sig.doli, sca.flav, chs.sord, …
## $ Magnetic.1 <fct> sca.rivu, sca.ghob, aca.duss, aca.bloc, …
## $ Whitsundays.1 <fct> sca.rivu, sig.doli, sca.flav, sca.nigr, …
## $ Keppels.1 <fct> sca.rivu, sca.ghob, sig.arge, sig.fusc, …
## $ X.7 <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
## $ X.8 <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
## $ X75..cumulative.minimum <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
## $ FG.1 <fct> Grazers, , , , , , , Planktivores, , , ,…
## $ Palms.2 <fct> sca.rivu, sig.doli, sca.flav, chs.sord, …
## $ Magnetic.2 <fct> sca.rivu, sca.ghob, aca.duss, aca.bloc, …
## $ Whitsundays.2 <fct> sca.rivu, sig.doli, sca.flav, sca.nigr, …
## $ Keppels.2 <fct> sca.rivu, sca.ghob, sig.arge, sig.fusc, …Exploratory data analyses
Total density
EDA_histograms(var='REGION', dat=fish, var.lookup=var.lookup)
NTR Polled
EDA_histograms(var='NTR.Pooled', dat=fish, var.lookup=var.lookup)
EDA_histograms(var='NTR.Pooled', group='REGION',dat=fish, var.lookup=var.lookup)
LHC % (live hard coral cover)
EDA_histograms(var='LHC', dat=fish, var.lookup=var.lookup)
EDA_density(var='LHC', group='REGION', dat=fish, var.lookup=var.lookup)
SC % (soft coral)
EDA_histograms(var='SC', dat=fish, var.lookup=var.lookup)
EDA_density(var='SC', group='REGION', dat=fish, var.lookup=var.lookup)
MA % (macroalgal cover)
EDA_histograms(var='MA', dat=fish, var.lookup=var.lookup)
EDA_density(var='MA', group='REGION', dat=fish, var.lookup=var.lookup)
Turf %
EDA_histograms(var='TURF', dat=fish, var.lookup=var.lookup)
EDA_density(var='TURF', group='REGION', dat=fish, var.lookup=var.lookup)
Unconsolidated %
EDA_histograms(var='UC', dat=fish, var.lookup=var.lookup)
EDA_density(var='UC', group='REGION', dat=fish, var.lookup=var.lookup)
Benthic richness
EDA_histograms(var='BR', dat=fish, var.lookup=var.lookup)
EDA_density(var='BR', group='REGION', dat=fish, var.lookup=var.lookup)
Coral morphological diversity
EDA_histograms(var='CMD', dat=fish, var.lookup=var.lookup)
EDA_density(var='CMD', group='REGION', dat=fish, var.lookup=var.lookup)
Slope
EDA_histograms(var='SLOPE', dat=fish, var.lookup=var.lookup)
EDA_density(var='SLOPE', group='REGION', dat=fish, var.lookup=var.lookup)
Rugosity
EDA_histograms(var='RUG', dat=fish, var.lookup=var.lookup)
EDA_density(var='RUG', group='REGION', dat=fish, var.lookup=var.lookup)
SCI (structural complexity)
EDA_histograms(var='SCI', dat=fish, var.lookup=var.lookup)
EDA_density(var='SCI', group='REGION', dat=fish, var.lookup=var.lookup)
Chlorophyl-a
EDA_histograms(var='CHL', dat=fish, var.lookup=var.lookup)
EDA_density(var='CHL', group='REGION', dat=fish, var.lookup=var.lookup)
Kd490
EDA_histograms(var='KD490', dat=fish, var.lookup=var.lookup)
EDA_density(var='KD490', group='REGION', dat=fish, var.lookup=var.lookup)
SST mean
EDA_histograms(var='SSTMEAN', dat=fish, var.lookup=var.lookup)
EDA_density(var='SSTMEAN', group='REGION', dat=fish, var.lookup=var.lookup)
SST anomoly
EDA_histograms(var='SSTANOM', dat=fish, var.lookup=var.lookup)
EDA_density(var='SSTANOM', group='REGION', dat=fish, var.lookup=var.lookup)
Wave exposure
EDA_histograms(var='WAVE', dat=fish, var.lookup=var.lookup)
EDA_density(var='WAVE', group='REGION', dat=fish, var.lookup=var.lookup)
Corrected Depth
EDA_histograms(var='DEPTH', dat=fish, var.lookup=var.lookup)
EDA_density(var='DEPTH', group='REGION', dat=fish, var.lookup=var.lookup)
Max DHW
EDA_histograms(var='DHW', dat=fish, var.lookup=var.lookup)
EDA_density(var='DHW', group='REGION', dat=fish, var.lookup=var.lookup)
Cyclone
EDA_histograms(var='CYCLONE', dat=fish, var.lookup=var.lookup)
EDA_density(var='CYCLONE', group='REGION', dat=fish, var.lookup=var.lookup)
Explore to primary water (weeks)
EDA_histograms(var='EXP', dat=fish, var.lookup=var.lookup)
EDA_density(var='EXP', group='REGION', dat=fish, var.lookup=var.lookup)
Prey density (for Plectropomus density and legal density)
EDA_histograms(var='PREY.DENSITY', dat=fish, var.lookup=var.lookup)
EDA_density(var='PREY.DENSITY', group='REGION', dat=fish, var.lookup=var.lookup)
Prey biomass (for Plectropomus density and legal density)
EDA_histograms(var='PREY.BIOMASS', dat=fish, var.lookup=var.lookup)
EDA_density(var='PREY.BIOMASS', group='REGION', dat=fish, var.lookup=var.lookup)
Results
formulas = list(
all = ~REGION + NTR.Pooled + LHC + SC + MA + TURF + UC + BR + CMD + SLOPE + RUG + SCI + CHL + KD490 + SSTANOM + WAVE + DEPTH + DHW + CYCLONE + EXP,
Palm = ~ NTR.Pooled + LHC + SC + MA + TURF + UC + BR + CMD + SLOPE + RUG + SCI + CHL + KD490 + SSTMEAN + SSTANOM + WAVE + DEPTH + DHW + CYCLONE + EXP,
Magnetic = ~ NTR.Pooled + LHC + SC + MA + TURF + UC + BR + CMD + SLOPE + RUG + SCI + CHL + KD490 + SSTMEAN + SSTANOM + WAVE + DEPTH + DHW + CYCLONE + EXP,
Whitsunday = ~ NTR.Pooled + LHC + SC + MA + TURF + UC + BR + CMD + SLOPE + RUG + SCI + CHL + KD490 + SSTMEAN + SSTANOM + WAVE + DEPTH + DHW + CYCLONE + EXP,
Keppel = ~ NTR.Pooled + LHC + SC + MA + TURF + UC + BR + CMD + SLOPE + RUG + SCI + CHL + KD490 + SSTMEAN + SSTANOM + WAVE + DEPTH + DHW + CYCLONE + EXP
)load('data/fish.RData')
load('data/var.lookup.RData')
resp.lookup = var.lookup %>% filter(Type=='Response') %>% droplevels
pred.lookup = var.lookup %>% filter(Type=='Predictor') %>% droplevels
groupings.all = vector('list', nrow(pred.lookup))
names(groupings.all) = pred.lookup$Abbreviation
groupings.region = vector('list', nrow(pred.lookup))
names(groupings.region) = pred.lookup$Abbreviation
for (i in 1:nrow(pred.lookup)) {
pred = as.character(pred.lookup[i,'Abbreviation'])
groupings.all[[pred]] = ifelse(pred=='REGION',NA,'REGION')
groupings.region[[pred]] = ifelse(pred=='NTR.Pooled',NA,'NTR.Pooled')
}
groupings.all = do.call('c',groupings.all)
groupings.region = do.call('c',groupings.region)
gfun = function(f, form) {
#print(f)
#print(form[[f]])
form = attr(terms(form[[f]]),'term.labels')
if (f=='all') return(groupings.all[form])
else return(groupings.region[form])
}
analyses = vector('list',nrow(resp.lookup))
names(analyses) = resp.lookup$Abbreviation
for (i in 1:nrow(resp.lookup)) {
resp = as.character(resp.lookup[i,'Abbreviation'])
fun = as.character(resp.lookup[i,'Transform'])
fam = as.character(resp.lookup[i,'Family'])
analyses[[resp]] = list('formulas' = lapply(formulas, function(f) update(f, paste0(fun,'(',resp,') ~.'))),
'family' = fam)
if (resp %in% c('PI','PTD','PLD')) analyses[[resp]][['formulas']] = lapply(analyses[[resp]][['formulas']], function(f) f=update(f, .~.+PREY.DENSITY))
if (resp %in% c('PTB','PLB')) analyses[[resp]][['formulas']] = lapply(analyses[[resp]][['formulas']], function(f) f=update(f, .~.+PREY.BIOMASS))
analyses[[resp]][['groups']] = sapply(c('all','Palm','Magnetic','Whitsunday','Keppel'), gfun, form=analyses[[resp]][['formulas']], USE.NAMES = TRUE,simplify=FALSE)
}for (a in 1:length(analyses)) {
resp=names(analyses)[a]
print(paste('Response =',resp))
for (f in 1:length(analyses[[a]]$formulas)) {
mod.name = names(analyses[[a]]$formulas)[f]
print(paste('Model =',mod.name))
MONOTONE = assignMonotone(fish, analyses[[a]]$formulas[[f]])
if (mod.name=='all') {fish.sub=fish
} else {fish.sub = fish %>% filter(REGION==mod.name)}
if (any(fish.sub[,as.character(get_response(analyses[[a]]$formulas[[f]]))]==0)) {
val = fish.sub[, as.character(get_response(analyses[[a]]$formulas[[f]]))]
val=min(val[val>0], na.rm=TRUE)
fish.sub[, as.character(get_response(analyses[[a]]$formulas[[f]]))] = fish.sub[,as.character(get_response(analyses[[a]]$formulas[[f]]))] + val
}
set.seed(123)
mod = abt(analyses[[a]]$formulas[[f]], data=fish.sub, distribution=analyses[[a]]$family,
cv.folds=10,interaction.depth=10,n.trees=10000, shrinkage=0.001, n.minobsinnode=2,
var.monotone=as.vector(MONOTONE))
p=plot.abts(mod, var.lookup, center=FALSE, type='response', return.grid=TRUE,
groupby=na.omit(unique(analyses[[a]]$groups[[f]])))#
thresholds=p$thresholds
save(thresholds, file=paste0('data/thresholds_',mod.name,'_',resp,'.RData'))
ps = p[['ps']]
p = p[['p']]
if (length(levels(fish$REGION))>1 || length(levels(fish$NTR.Pooled))>1) {
p = common_legend(p)
ps = common_legend(ps)
}
## version for the supplimentary
#do.call('grid.arrange', p) ## version for the supplimentary
ggsave(filename=paste0('output/figures/data.all.abt.',mod.name,'_',resp,'_ABT.png'), do.call('grid.arrange', p), width=15, height=10, dpi=300)
ggsave(filename=paste0('output/figures/data.all.abt.',mod.name,'_',resp,'_ABT.pdf'), do.call('grid.arrange', p), width=15, height=10)
pred.1=stats.abt(mod, fitMethod=1, analysis = analyses[[a]]$groups[[f]])
save(pred.1, file=paste0('data/pred.1_',mod.name,'_',resp,'.RData'))
optim=summarize_values(pred.1$optim, type='optim') %>%
dplyr::rename_at(vars(-contains('Var')), paste0, '.optim')
rel.inf = summarize_values(pred.1$rel.imp, type='Rel.inf') %>%
dplyr::rename_at(vars(-contains('Var')), paste0, '.rel.inf')
R2=summarize_values(pred.1$R2.value, 'R2') %>%
dplyr::rename_at(vars(-contains('Var')), paste0, '.R2')
stats = optim %>% full_join(R2) %>% full_join(rel.inf) %>%
left_join(var.lookup %>% dplyr::select(Var=Abbreviation, pretty.name))
save(stats, file=paste0('data/stats_',mod.name,'_',resp,'.RData'))
write.csv(stats %>% as.data.frame, file=paste0('output/data/stats_',mod.name,'_',resp,'.csv'), quote=FALSE, row.names=FALSE)
## Version full model with just the relative importance and the substantial panels
ps1 = arrangeGrob(ps[[1]],ps[[length(ps)]], widths=c(3,1))
numberOfPlots = length(ps)-2 #minus 1 since one of the items is the legend and minus one for the relative importance plot
numberOfPlotRows = ceiling(numberOfPlots / 2)
g=arrangeGrob(ps1, do.call('arrangeGrob', c(ps[c(-1, -length(ps))], list(ncol=2))), heights=c(2,numberOfPlotRows))
ggsave(filename=paste0('output/figures/data.all.abt.',mod.name,'_',resp,'_ABT_short.png'), g, width=10, height=2+(1.7*numberOfPlotRows), dpi=300)
ggsave(filename=paste0('output/figures/data.all.abt.',mod.name,'_',resp,'_ABT_short.pdf'), g, width=10, height=2+(1.7*numberOfPlotRows))
## Version with a grid of partials in the lower right corner of the influence figure
## The location and size of the figure in figure will be determined by the scale of relative influence
## ideally, we want the subfigure to be 3/4 of the width and 3/4 of the height
ymax=max(pretty(as.numeric(as.character(stats$upper.rel.inf))))
ymin=ymax*1/4
xmax=(length(p)-2)*3/4
if (ymin <= 100/(length(p)-2)) { # if the left side of the subfigure is too close to the dashed vertical line
ymax = ymax + 5
ymin=ymax*1/4
ps[[1]] = ps[[1]] + scale_y_continuous('Relative Importance', limits=c(0,ymax))
}
ps[2:(length(ps)-1)]=lapply(ps[2:(length(ps)-1)], function(f) f+theme(axis.title.y=element_blank()))
ps2=do.call('arrangeGrob', c(ps[c(length(ps),2:(length(ps)-1))], list(ncol=2)))
#g= ps[[1]] + annotation_custom(grob=ps2, xmin=1, xmax=15, ymin=10, ymax=Inf)
g= ps[[1]] + annotation_custom(grob=ps2, xmin=1, xmax=xmax, ymin=ymin, ymax=Inf)
ggsave(filename=paste0('output/figures/data.all.abt.',mod.name,'_',resp,'_ABT_short_in.png'), g, width=10, height=8, dpi=300)
ggsave(filename=paste0('output/figures/data.all.abt.',mod.name,'_',resp,'_ABT_short_in.pdf'), g, width=10, height=8)
## Version for vertical table of regional relative importance and substantial panels
g=do.call('arrangeGrob', c(ps[1:4], list(ncol=1, heights=c(1.5, rep(1,3)))))
#save(g, file=paste0('data/g.abt_',mod_num,'_',resp,'.RData'))
#save(data.all.abt, file=paste0('data/data.all.abt_',mod_num,'_',resp,'.RData'))
#save(rel.importance, file=paste0('data/rel.importance_',mod_num,'_',resp,'.RData'))
#save(thresholds, file=paste0('data/thresholds_', mod_num,'_',resp,'.RData'))
#write.table(bind_rows(thresholds, .id='Var'), file=paste0('data/thresholds_',mod_num,'_',resp,'.csv'), quote=FALSE, row.names=FALSE)
rm(mod,pred.1,p,ps,ps1,g,stats,optim,rel.inf,R2,fit)
gc()
}
}
#summary(mod)
#plot(mod,c(15,1))load('data/var.lookup.RData')
resp.lookup = var.lookup %>% filter(Type=='Response') %>% droplevels
mods = names(formulas)
for (a in 1:length(analyses)) {
resp=names(analyses)[a]
cat(paste('## ',resp.lookup$pretty.name[resp.lookup$Abbreviation==resp],' {.tabset .tabset-pills} \n\n'))
for (f in 1:length(analyses[[a]]$formulas)) {
mod.name = names(analyses[[a]]$formulas)[f]
cat(paste('### ',mod.name,'\n\n'))
cat(paste0('\n\n'))
## Now for a tables
load(file=paste0('data/stats_',mod.name,'_',resp,'.RData'))
pretty.stats(stats, var.lookup) %>% rmarkdown:::paged_table_html() %>% cat
load(file=paste0('data/thresholds_',mod.name,'_',resp,'.RData'))
pretty.thresholds(thresholds, stats, var.lookup) %>% rmarkdown:::paged_table_html() %>% cat
}
}Total density
all
Palm
Magnetic
Whitsunday
Keppel
Total species richness
all
Palm
Magnetic
Whitsunday
Keppel
Benthic invertivores
all
Palm
Magnetic
Whitsunday
Keppel
Grazers
all
Palm
Magnetic
Whitsunday
Keppel
Corallivores
all
Palm
Magnetic
Whitsunday
Keppel
Omnivores
all
Palm
Magnetic
Whitsunday
Keppel
Planktivores
all
Palm
Magnetic
Whitsunday
Keppel
Carnivores
all
Palm
Magnetic
Whitsunday
Keppel
Piscivores
all
Palm
Magnetic
Whitsunday
Keppel
Farmers
all
Palm
Magnetic
Whitsunday
Keppel
Plectropomus total density
all
Palm
Magnetic
Whitsunday
Keppel
Plectropomus total biomass
all
Palm
Magnetic
Whitsunday
Keppel
Plectropomus legal density
all
Palm
Magnetic
Whitsunday
Keppel
Plectropomus legal biomass
all
Palm
Magnetic
Whitsunday
Keppel
